# Data Manipulation and Linear Algebra
import pandas as pd
import numpy as np

# Machine Learning and Data Preprocessing
from sklearn.decomposition import PCA, TruncatedSVD, NMF
from sklearn.preprocessing import StandardScaler

# Data Visualization
import matplotlib.pyplot as plt
import seaborn as sns

# Web Interaction and Display
from IPython.display import Image, display, HTML

# Warning Handling
import warnings
from sklearn.exceptions import ConvergenceWarning
warnings.filterwarnings("ignore", category=ConvergenceWarning)

# Additional JavaScript for toggling code display in Jupyter Notebooks
HTML(
    """
<script src="https://cdnjs.cloudflare.com/ajax/libs/jquery/2.0.3/jquery.min.js "></script>
<script>
code_show=true; 
function code_toggle() {
 if (code_show){
 $('div.jp-CodeCell > div.jp-Cell-inputWrapper').hide();
 } else {
 $('div.jp-CodeCell > div.jp-Cell-inputWrapper').show();
 }
 code_show = !code_show
} 
$( document ).ready(code_toggle);
</script>
<form action="javascript:code_toggle()"><input type="submit"
    value="Click here to toggle on/off the raw code."></form>
"""
)

# Artist Count in top tracks
df_tracks = pd.read_csv("top_tracks.csv")
artist = df_tracks["artist"].value_counts()[:5].to_frame()
plt.figure(figsize=(10, 6), dpi=250)
plt.barh(
    artist.index,
    artist["count"],
    color="#1ccf54",
)
plt.xlabel("Count")
plt.xticks(rotation=0, ha="right")
plt.title("Artist Count")
plt.show()

# Artist Count in the track pool
df_pool = pd.read_csv("track_pool.csv")
artist = df_pool["artist"].value_counts()[:5].to_frame()
plt.figure(figsize=(10, 6), dpi=250)
plt.barh(artist.index, artist["count"], color="#1ccf54")
plt.xlabel("Count")
plt.xticks(rotation=0, ha="right")
plt.title("Artist Count")
plt.show()

# Correl Plot of Audio Features
num = df_pool.select_dtypes(exclude="object").columns
num_df = df_pool[num]
correlation_matrix = num_df.corr()
plt.figure(figsize=(10, 8), dpi=250)
sns.heatmap(correlation_matrix, cmap="YlGnBu", linewidths=0.5)
plt.title("Correlation Heatmap of Audio Features")
plt.show()

audio_features = [
    "danceability",
    "energy",
    "loudness",
    "speechiness",
    "acousticness",
    "instrumentalness",
    "liveness",
    "valence",
    "tempo",
]
df_pool = df_pool[audio_features]

# Create a single figure with 9 subplots (one for each audio feature)
plt.figure(figsize=(18, 12))
for i, feature in enumerate(audio_features, start=1):
    plt.subplot(3, 3, i)
    sns.boxplot(data=df_pool[feature], color="#1ccf54", width=0.5,
                linewidth=1.5)
    plt.title(feature.title(), fontsize=14)

plt.tight_layout()
plt.show()

def clean_audio_features(consolidated_audio_features):
    """
    Cleans and preprocesses audio feature data.

    Parameters
    ----------
    consolidated_audio_features : DataFrame
        A pandas DataFrame containing the audio features.

    Returns
    -------
    DataFrame
        A cleaned DataFrame with normalized audio features, excluding 'key'
        and 'mode' features.
    """
    # Retain only the audio features without 'mode' and 'key'
    unwanted_features = ["key", "mode"]
    wanted_features = [
        feature
        for feature in consolidated_audio_features
        if feature not in unwanted_features
    ]
    cleandf_audio_features = consolidated_audio_features.reindex(
        columns=wanted_features
    )

    # Get only the columns of audio features
    cleandf_audio_features = cleandf_audio_features.iloc[:, :9]
    wanted_features = list(cleandf_audio_features.columns)

    # All matrix entries must be nonnegative
    cleandf_audio_features = np.abs(cleandf_audio_features)

    # Scale the features to avoid over/underemphasizing of features
    # Use min-max scaling
    from sklearn.preprocessing import MinMaxScaler

    min_max_scaler = MinMaxScaler()
    cleandf_audio_features = min_max_scaler.fit_transform(
        cleandf_audio_features)
    cleandf_audio_features = pd.DataFrame(
        cleandf_audio_features, columns=wanted_features
    )

    return cleandf_audio_features

# Track Pool for NMF 
consolidated_audio_features = pd.read_csv('track_pool.csv')
df_audio_features = clean_audio_features(consolidated_audio_features)
df_audio_features.head()

# Initiating the NMF model
nmf_model = NMF(max_iter=10000)
A = np.array(df_audio_features)
W = nmf_model.fit_transform(A)
# W.shape Should be (number of archetypes, number of tracks)
H = nmf_model.components_.T
# H.shape Should be (number of archetypes, number of audio features)

for archetype, coefficients in enumerate(H):
    # 'coefficients' is the importance of each audio feature to an archetype
    # Using the coefficients, get the index of the top five audio features
    top_feature_indices = coefficients.argsort()[::-1][:5]

    # Get the corresponding audio feature of each index
    top_feature_names = df_audio_features.columns[top_feature_indices]

    print(f"Archetype {archetype+1} Top Features:")
    print(list(top_feature_names), end="\n\n")

# Add back the id of each track
audio_matrix = pd.concat(
    [df_audio_features.iloc[:, :11], consolidated_audio_features["id"]], axis=1
)

# Assign an archetype to each user using W
# Matrix entries of W refer to how much of each archetype is found in each track
# The highest matrix entry gives the archetype that is dominant in each track
audio_matrix["Archetype"] = (
    np.argmax(W, axis=1) + 1
)  # Add 1 to make the archetype labels start from 1
audio_matrix

Archetype 1 Top Features:
['valence', 'instrumentalness', 'energy', 'tempo', 'liveness']

Archetype 2 Top Features:
['liveness', 'danceability', 'instrumentalness', 'tempo', 'valence']

Archetype 3 Top Features:
['tempo', 'loudness', 'valence', 'liveness', 'instrumentalness']

Archetype 4 Top Features:
['danceability', 'tempo', 'loudness', 'valence', 'liveness']

Archetype 5 Top Features:
['energy', 'tempo', 'valence', 'liveness', 'instrumentalness']

Archetype 6 Top Features:
['loudness', 'danceability', 'energy', 'instrumentalness', 'tempo']

Archetype 7 Top Features:
['acousticness', 'tempo', 'energy', 'valence', 'liveness']

Archetype 8 Top Features:
['instrumentalness', 'speechiness', 'acousticness', 'tempo', 'valence']

Archetype 9 Top Features:
['speechiness', 'liveness', 'valence', 'tempo', 'instrumentalness']

# Sparsity Plot of NMF for Non-Zero Elements
audio_features = audio_matrix.columns.tolist()[:9]
fig, ax = plt.subplots(dpi=250)
ax.spy(H)
ax.set_xticks(range(len(audio_features)))
ax.set_yticks(range(len(audio_features)))
ax.set_yticklabels(df_audio_features.columns);

# Reconstruction Error Plot for NMF
n_components = list(range(1, 10))
errors = []
n_components_list = []
for n_component in n_components:
    nmf_ = NMF(n_component, max_iter=10000)
    nmf_.fit(A)
    n_components_list.append(n_component)
    errors.append(nmf_.reconstruction_err_)
plt.figure(dpi=250)
plt.plot(n_components, errors, "-o", color="#1ccf54")
plt.xlabel(r"$n_\text{components}$")
plt.ylabel("reconstruction error");

# Print each n_component and its corresponding reconstruction error
for n_components, error in zip(n_components_list, errors):
    print(f"Archetypes: {n_components}")
    print(f"Reconstruction error: {error:.6f}")
    print("-" * 31)

Archetypes: 1
Reconstruction error: 23.087236
-------------------------------
Archetypes: 2
Reconstruction error: 18.141172
-------------------------------
Archetypes: 3
Reconstruction error: 14.174652
-------------------------------
Archetypes: 4
Reconstruction error: 11.223809
-------------------------------
Archetypes: 5
Reconstruction error: 9.082671
-------------------------------
Archetypes: 6
Reconstruction error: 7.077012
-------------------------------
Archetypes: 7
Reconstruction error: 4.829430
-------------------------------
Archetypes: 8
Reconstruction error: 3.618093
-------------------------------
Archetypes: 9
Reconstruction error: 0.259035
-------------------------------

# Using wanted features for tracks
features = [
    "danceability",
    "energy",
    "loudness",
    "speechiness",
    "acousticness",
    "instrumentalness",
    "liveness",
    "valence",
    "tempo",
]
df_features = df_pool[features]

# Scaling data
scaler = StandardScaler()
df_scaled = scaler.fit_transform(df_features)

# Fitting on all components
pca = PCA()
X_new = pca.fit_transform(df_scaled)
variance_explained = pca.explained_variance_ratio_

def de_correlate_df(df):
    """
    De-correlate the columns of a DataFrame by shuffling each column
    independently.

    Parameters
    ----------
    df: pandas.DataFrame
        The input DataFrame to de-correlate.

    Returns
    -------
    X_aux: pandas.DataFrame
        A new DataFrame with the same columns as the input DataFrame,
        but with each column's values shuffled independently.

    """
    X_aux = df.copy()
    for col in df.columns:
        X_aux[col] = df[col].sample(len(df)).values
    return X_aux

# Define the number of tests and create a blank matrix
N_permutations = 1000
variance = np.zeros((N_permutations, len(features)))
pca_iter = PCA()
# Generate the new datasets and save the result
for i in range(N_permutations):
    X_aux = de_correlate_df(pd.DataFrame(df_scaled))
    pca_iter.fit(X_aux)
    variance[i, :] = pca_iter.explained_variance_ratio_
p_val = np.sum(variance > variance_explained, axis=0) / N_permutations

# Plot of P-values of Permutation Test
plt.figure(figsize=(4, 2), dpi=250)
plt.plot(range(1, len(variance_explained) + 1), p_val, "o-", color="#1ccf54")
plt.ylim(-0.1, 1.1)
plt.xlabel("Principal Components")
plt.ylabel("p-value")
plt.show()

# Plot for cumulative and single Variance Explained
plt.figure(dpi=250)
plt.plot(
    range(1, len(variance_explained) + 1), variance_explained, "o-", color="#1ccf54"
)
plt.plot(
    range(1, len(variance_explained) + 1),
    variance_explained.cumsum(),
    "o-",
    color="#f465bf",
)
plt.ylim(0, 1)
plt.xlabel("PC")
plt.ylabel("variance explained")
plt.legend(["Single", "Cumulative"], loc=0);

# Extracting the Principal Components
n_components = pca.components_[:3].T
pc = pd.DataFrame(
    n_components,
    index=features,
    columns=["PC" + str(x) for x in range(1, n_components.shape[1] + 1)],
)

# Plotting coefficients of features projected to PC1
colors = ["gray"] * len(features)
colors[4] = "#65D46E"
colors[6] = "#65D46E"
colors[7] = "#65D46E"
plt.figure(figsize=(10, 6), dpi=250)
pc["PC1"][::-1].plot(kind="barh", color=colors)
plt.title("PRINCIPAL COMPONENT 1\n\n" + "←INTENSE vs MELLOW→")
plt.xlabel("Coefficients")
plt.ylabel("Features")
plt.xlim(-0.6, 0.6)
plt.show()

# Plotting coefficients of features projected to PC2
colors = ["gray"] * len(features)
colors[1] = "#65D46E"
colors[2] = "#65D46E"
colors[8] = "#65D46E"
plt.figure(figsize=(10, 6), dpi=250)
pc["PC2"][::-1].plot(kind="barh", color=colors)
plt.title("PRINCIPAL COMPONENT 2\n\n" + "←HIGH-SPIRITED vs MELANCHOLIC→")
plt.xlabel("Coefficients")
plt.ylabel("Features")
plt.xlim(-0.6, 0.6)
plt.show()

# Plotting coefficients of features projected to PC3
colors = ["gray"] * len(features)
colors[2] = "#65D46E"
colors[3] = "#65D46E"
colors[5] = "#65D46E"
plt.figure(figsize=(10, 6), dpi=250)
pc["PC3"][::-1].plot(kind="barh", color=colors)
plt.title("PRINCIPAL COMPONENT 3\n\n" + "←NON-LYRICAL vs LYRICAL→")
plt.xlabel("Coefficients")
plt.ylabel("Features")
plt.xlim(-0.6, 0.6)
plt.show()

# Fitting on optimal number of components
pca = PCA(3)
X_new = pca.fit(df_scaled)
coord = pca.components_.T

# Fitting User to optimal PCA components 
df_test = pd.read_csv("top_tracks.csv")
df_test_scaled = scaler.transform(df_test[features])
result = pca.transform(df_test_scaled)
df_result = pd.DataFrame(result, columns=["PC1", "PC2", "PC3"])
df_result["Archetype"] = (
    df_result["PC1"].apply(lambda x: 1 if x > 0 else 0)
    + df_result["PC2"].apply(lambda x: 3 if x > 0 else 0)
    + df_result["PC3"].apply(lambda x: 5 if x > 0 else 0)
)
mapper = {
    9: "Poignant Poet",
    4: "Melancholic Melodist",
    6: "Ambient Alchemist",
    1: "Serenity Savant",
    8: "Vocal Voyager",
    3: "Reflective Rythmist",
    5: "Lyric Luminary",
    0: "High-Spirited Harmonist",
}
df_result["Archetype"] = df_result["Archetype"].map(mapper)
df_result["user"] = df_test["user"].to_list()
most_frequent = (df_result.groupby("user")["Archetype"]).agg(
    lambda x: x.mode().iloc[0])
user_index = [f"User {i}" for i in range(1, 6)]
most_frequent.index = user_index

most_frequent

User 1          Ambient Alchemist
User 2          Ambient Alchemist
User 3             Lyric Luminary
User 4    High-Spirited Harmonist
User 5          Ambient Alchemist
Name: Archetype, dtype: object

Attribute	Type	Description
acousticness	float	A confidence measure from 0.0 to 1.0 of whether the track is acoustic. 1.0 represents high confidence the track is acoustic.
danceability	float	Danceability describes how suitable a track is for dancing based on a combination of musical elements including tempo, rhythm stability, beat strength, and overall regularity. A value of 0.0 is least danceable and 1.0 is most danceable.
duration_ms	integer	The duration of the track in milliseconds.
energy	float	Energy is a measure from 0.0 to 1.0 and represents a perceptual measure of intensity and activity. Typically, energetic tracks feel fast, loud, and noisy.
instrumentalness	float	The closer the instrumentalness value is to 1.0, the greater likelihood the track contains no vocal content. Values above 0.5 are intended to represent instrumental tracks, but confidence is higher as the value approaches 1.0.
key	integer	The converted mana cost of the card.
liveness	float	Detects the presence of an audience in the recording. Higher liveness values represent an increased probability that the track was performed live. A value above 0.8 provides strong likelihood that the track is live.
loudness	float	The overall loudness of a track in decibels (dB). Loudness values are averaged across the entire track and are useful for comparing relative loudness of tracks.
mode	integer	Mode indicates the modality (major or minor) of a track, the type of scale from which its melodic content is derived. Major is represented by 1 and minor is 0.
speechiness	float	Speechiness detects the presence of spoken words in a track. The more exclusively speech-like the recording (e.g. talk show, audio book, poetry), the closer to 1.0 the attribute value.
tempo	float	The overall estimated tempo of a track in beats per minute (BPM). In musical terminology, tempo is the speed or pace of a given piece and derives directly from the average beat duration.
time_signature	integer	An estimated time signature. The time signature (meter) is a notational convention to specify how many beats are in each bar (or measure). The time signature ranges from 3 to 7 indicating time signatures of "3/4", to "7/4".
valence	float	A measure from 0.0 to 1.0 describing the musical positiveness conveyed by a track. Tracks with high valence sound more positive (e.g. happy, cheerful, euphoric), while tracks with low valence sound more negative (e.g. sad, depressed, angry).

	danceability	energy	loudness	speechiness	acousticness	instrumentalness	liveness	valence	tempo
0	0.748892	0.807101	0.249543	0.015001	0.044778	0.000041	0.034801	0.793684	0.362562
1	0.559141	0.275624	0.283591	0.008526	0.761044	0.003582	0.079663	0.287368	0.574333
2	0.521645	0.331887	0.291086	0.038420	0.157630	0.930612	0.050919	0.596842	0.305260
3	0.596637	0.206301	0.218620	0.008526	0.714859	0.000000	0.127400	0.295789	0.514339
4	0.676173	0.189222	0.323254	0.027736	0.949799	0.000007	0.331691	0.431579	0.483237

	danceability	energy	loudness	speechiness	acousticness	instrumentalness	liveness	valence	tempo	id	Archetype
0	0.748892	0.807101	0.249543	0.015001	0.044778	0.000041	0.034801	0.793684	0.362562	0obzwoGjXdORtAwYWZ9s53	6
1	0.559141	0.275624	0.283591	0.008526	0.761044	0.003582	0.079663	0.287368	0.574333	4oBEeMxkuW4t7a8qIupetG	8
2	0.521645	0.331887	0.291086	0.038420	0.157630	0.930612	0.050919	0.596842	0.305260	1BMep2eRJLHOZnpL8Kd0lY	6
3	0.596637	0.206301	0.218620	0.008526	0.714859	0.000000	0.127400	0.295789	0.514339	7noCkklUhEAoj8GJkbAuHq	8
4	0.676173	0.189222	0.323254	0.027736	0.949799	0.000007	0.331691	0.431579	0.483237	1SZ63wN0pk18Dr4Epyhcsf	8
...	...	...	...	...	...	...	...	...	...	...	...
1246	0.539825	0.172142	0.462217	0.019426	0.824297	0.000110	0.102762	0.498947	0.307914	3Fn6K2roSmzhKzp5BBJKJl	9
1247	0.721623	0.616212	0.354685	0.033887	0.653614	0.986735	0.093522	0.662105	0.398669	2G3CRca8Si7vrCeDain6d0	6
1248	0.358027	0.542870	0.205763	0.012519	0.283132	0.505102	0.093522	0.283158	0.629315	6Ru1MJhTTFFnbSUqTGmrtM	7
1249	0.738666	0.683525	0.297718	0.030218	0.038352	0.000337	0.014783	0.798947	0.637420	2yrAD69v6rP6bZhIjs7kOy	6
1250	0.602318	0.371069	0.279398	0.015325	0.790160	0.005449	0.064983	0.722105	0.449578	2ARHr1CgtVWdRhlOIuBD61	6

ABSTRACT

INTRODUCTION

Background

Problem Statement

Objectives

DATA SOURCES AND DESCRIPTION

METHODOLOGY

EXPLORATORY DATA ANALYSIS

DIMENSIONALITY REDUCTION

NMF¶

PCA¶

Permutation Test¶

RESULTS AND DISCUSSION

DECIDING ARCHETYPES

Interpretation of Principal Components¶

PREDICTION

CONCLUSION

RECOMMENDATION

References¶

Appendix¶

Description of the Archetypes¶